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I. INTRODUCTION 

Relaxation, the return of a perturbed system into equi- 
librium, is one of the most fundamental processes of phys- 
ical systems. The simplest relaxation model is the expo- 
nential or Debye law exp(— i/r) incorporating the relax- 
ation time scale r. In many cases, however, significant 
deviations from the exponential law have been observed. 
The two most prominent generalized relaxation laws to 
accommodate these deviations are the stretched expo- 
nential law exp(— [t/r] a ) with < a < 1 and the asymp- 
totic power or Pareto law (sometimes also called Nutting 
law) (1 + [i/ T F) ( ~ 1 " 1 witn /3 > 0. When compared to the 
simple exponential law, both these generalized laws cor- 
respond to a broad distribution of relaxation times 
Generalized relaxation patterns are observed on many 
scales, ranging from single molecules 043 over dielectric 
response @ to macroscopic viscoelasticity 

Theoretical approaches to describe generalized re- 
laxation processes include parallel relaxation channels 
(ToL [TlT | , hierarchically constrained dynamics giving rise 
to complex serial relaxation [ljj], and defect diffusion 
models [l3j]. For the stretched exponential law Klaftcr 
and Shlesinger demonstrated the common universal fea- 
tures behind these approaches [14j , There exist also ex- 
tensions of stretched exponentials in models of dynamic 
relaxation channels Yet another approach to gen- 

eralized relaxation dynamics is that of fractional-order, 
viscoelastic mechanical bodies [l6l - l2fjj ; for their histor- 
ical development see Ref. [2l[ . Similar methods based on 
generalized dynamic equations and stochastic approaches 
based on stable distributions have been used to describe 
dielectric relaxation behavior 122, ^3 . 



In this paper we establish a robust Poissonian ap- 
proach to complex relaxation processes based on indi- 
vidual reactions between excitations, that are randomly 
scattered around a reaction center in a general embed- 
ding space. The resulting model is termed RARE, the ac- 
ronym standing for RAndom RElaxations. The model's 
input quantities are the spatial scattering statistics of the 
excitations around the reaction center, and the chemical 
reaction rates between the excitations and the reaction 
center as a function of their mutual distance. The RARE 
model is considerably general, it has a robust framework, 
and it allows for a fairly intuitive interpretation of com- 
plex relaxation processes. The RARE model is analyzed 
in detail, and analytic results regarding the duration and 
the range of the model's random relaxation processes are 
obtained in closed form. 

The paper is organized as follows. The RARE model 
is introduced, intuitively explained, and rigorously con- 
structed in section II. The detailed stochastic analysis of 
the RARE model and a Monte-Carlo algorithm for the 
simulation of the model's random relaxation processes 
are presented, respectively, in sections III and IV. The 
results presented in section III establish a comprehensive 
statistical picture of the duration and the range of the 
model's random relaxation processes: marginal distribu- 
tions, joint distribution, and conditional distributions. 
The case of power-law inputs, which are shown to yield 
stretched exponential relaxation patterns and asymptot- 
ically Paretian relaxation ranges, is addressed in section 
V. The thermodynamic limit of the RARE model is in- 
vestigated in section VI. Detailed proofs of the results 
stated along the paper are given in the Appendices. 
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The RARE model we introduce and explore in this pa- 
per is described as follows. A reaction center is placed 
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at an arbitrary point P of a general metric space Ai, 
and a countable collection of excitations is scattered ran- 
domly across the space M.. The excitations are labeled 
with the index i, and the position of excitation i is the 
random point Pj. The distance in the metric space M. 
is measured by a general metric function d(-, •), and the 
distance between the reaction center and excitation i is 
Di = d(Pi,P). Excitation i is equipped with a random 
timer Ti , and a reaction between the center and the excit- 
ations occurs upon the first timer-expiration event. The 
model is statistically characterized by two random vari- 
ables, the reaction time T and the reaction range X, 
which arc defined as follows: The reaction time T is the 
time elapsing until the first timer expires: 



T = min{Ti}. 



(1) 



The reaction range X is the distance between the reac- 
tion center and the excitation whose timer first expired: 



X = ^AI(T = Ti). 



(2) 



In Eq.[2]and hereinafter, 1(E) denotes the indicator func- 
tion of an event E (i.e., 1(E) = 1 if the event did occur, 
and 1(E) = if the event did not occur). 

Setting the space M. in which the reaction takes place 
to be a metric space yields great versatility. Indeed, M 
can be an Euclidean space of arbitrary dimension, a non- 
Euclidean space such as an elliptic space or a hyperbolic 
space, a general surface or landscape, a fractal object, a 
network, etc. We now turn to specify the RARE model 
assumptions, and thereafter present a preliminary ana- 
lysis. 



A. Model assumptions 

To quantify the RARE model the distributions of the 
random points {Pi}, as well as the distribution of the 
random timers {Ti}, need to be specified. A highly ap- 
plicable and effective statistical methodology to model 
the random scatte ring of points in general spaces are 
Poisson processes [24j. Poisson processes have a wide 
span of applications ranging from insurance and finance 
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to que uing systems [26|, and from fractals 
power- laws [30]. We henceforth assume that the random 
points {Pi} form a general Poisson process defined on 
the metric space M.. Consequently, the 'displacement 
theorem' of the theory of Poisson processes (see section 
5.5 in Ref. 24]) implies that the distances {Di} form a 
general Poisson process defined on the positive half-line 
(0, oo). In what follows we denote by p(x) the average 
number of excitations, which are within a distance x of 
the reaction center, 



p(x) = E 



j2m<x) 



In Eq. Q and in the following, E denotes the math- 
ematical expectation. Namely, if £ is a real-valued ran- 
dom variable governed by the probability density func- 
tion f^(x) (x real), and (j>(x) is a real- valued function 
defined on the real line, then E[0(£)] = <p(x)f^(x)dx. 
Moreover, in Eq. Q the function p(x) is monotone in- 
creasing, and is assumed to start at zero: p(0) = 0. We 
note that the derivative p'(x) is the Poissonian intensity 
of the Poisson process {Di} [24|. 

As noted in the introduction, the most common stat- 
istical law for reactions in the physical sciences is the ex- 
ponential law. Given the distances {Di} we now assume 
that the timers {Ti} are independent and exponentially 
distributed random variables, and that the exponential 
distribution of timer Ti is determined by the distance 
Di. In what follows we denote by r)(x) the exponential 
rate of the timers as a function of the distance variable 
x. Namely, given the distance Di, the timer Ti is expo- 
nentially distributed with tail distribution function 



Pr(T, ; > t\Di) = exp{-T](Di)t}, t > 0. 



(4) 



x > 0. 



(3) 



In other words, given the distance Di, the timer Ti is 
exponentially distributed with mean E[Tj | Di) — l/r/(Di). 
Typically, the function rj(x) is monotone decreasing in 
the distance variable x. 

The RARE model is quantified by the pair of func- 
tions introduced above, the scattering function p(x) and 
the reactivity function rj(x). The scattering function p(x) 
quantifies the underlying spatial scattering of the excit- 
ations, and the reactivity function r/(x) quantifies the 
underlying distance-dependent reaction rate. The 'in- 
puts' of the RARE model are the scattering function p(x) 
and the reactivity function rj(x), and the 'outputs' of the 
RARE model are the reaction time T and the reaction 
range X. In what follows we establish the statistics of the 
random outputs based on the given deterministic inputs. 

The RARE model introduced herein can be viewed as 
a descendant of a unified donor-acceptor energy transfer 
model presented by Blumen [28]. In its basic form this 
model considers an iV-site lattice, in which each site is oc- 
cupied by an excitation with probability p (independent 
of all other sites), and the reactions between the excit- 
ations and the reaction center are governed by Eq. (TJJ. 
Thomas and coworkers (2(| consider a model similar to 
the one developed herein, but with a uniform distribution 
of excitations. The RARE model generalizes these mod- 
els, as it (i) implements the notion of Poisson processes, 
(ii) regards both the reaction time T and the reaction 
range X, and (iii) considers a non- uniform distribution 
of excitations. On the one hand, the implementation 
of a Poisson scattering of excitations in a general met- 
ric space allows for a high variability and robustness; on 
the other hand, this Poisson modeling is highly tract- 
able, as it quantifies all the spatial and scattering de- 
tails (both potentially having infinitely many degrees of 
freedom) into one single function, the scattering function 
p(x). Moreover, in this paper we provide a stochastic 
analysis of the reaction pair (T,X), whereas in Ref. [28[ 
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only the reaction time T was analyzed. 



B. Preliminary analysis 

In the next section we present a detailed stochastic 
analysis of the RARE model. To facilitate the stochastic 
analysis of section [TTT] we present here a preliminary ana- 
lysis. In what follows we set R(x) to denote the aggregate 
exponential rate corresponding to excitations, whose dis- 
tance from the reaction center is greater than x, 



R(x) = J2v(Di)I{Di >x), x>0. 



(5) 



Note that the aggregate rate R(x) is a stochastic process 
parameterized by the distance variable x. 

Conditioned on the realizations of distances {-Dj}, the 
aggregate rate R(x) yields compact formulas for the con- 
ditional distributions of the reaction time T and the re- 
action range X. Indeed, the assumptions of the RARE 
model, combined with the statistical properties of the 
minima of independent exponential random variables 
31], implies that: (i) The conditional distribution of the 
reaction time T is given by the tail distribution function 



Pr(T > t\{Di}) = exp{-fcR(0)}, t > 0. 



(6) 



(ii) The conditional distribution of the reaction range X 
is given by the tail distribution function 



Pr(A > x\{Di}) 



R(x) 

W)' 



x>0. 



(7) 



(iii) The joint conditional distribution of the reaction pair 
(T, X) is given by the joint tail distribution function 

Pr(T > t, X > x\{Di}) = exp{-i#(0)}S ; t ,x > 0. 

(8) 

Note that given the realizations of the distances {Di} 
the reaction time T and the reaction range X turn out to 
be independent random variables. Indeed, the joint tail 
distribution function of the reaction pair (T, X) equals 
the product of the tail distribution functions of the reac- 
tion time T and the reaction range X, 

Pr(T>t,X>x\{Di}) = 

Pr(T > *|{A})Pr(X > x\{Di}), t,x>0. (9) 



III. STOCHASTIC ANALYSIS 

The preliminary analysis presented in section |TT] 
provides us with the conditional distribution of the re- 
action pair (T, X), conditioned on the realizations of dis- 
tances {Di}. A stochastic analysis detailed the Appen- 
dices A-D shifts us from the aforementioned conditional 



distribution to the (unconditional) distribution of the re- 
action pair (T, X). In this section we present the key 
results of this stochastic analysis. 

Reaction time. The tail distribution function of the 
reaction time T is given by 



Pr(T >t) = exp 



[1 - exp{-r)(x)t}]p(dx] 



(10) 

for t > 0. In turn, taking the limit t oo in Eq. (|T0|) 
yields the probability that a reaction never takes place, 



Pr(T = oo) = exp{— p(oo)}. 



(11) 



Eq. (fTTj) implies that a reaction takes place with cer- 
tainty, i.e., Pr(T < oc) = 1, if and only if the scattering 
function p(x) diverges, which in turn occurs if and only if 
there are infinitely many excitations. Moreover, differen- 
tiating Eq. (fT0| with respect to the time variable t yields 
the probability density function of the reaction time T, 

poo 

f T (t) = Pr(T >t) exp{-rj(x)t}r)(x)p(dx), t > 0. 
Jo 

(12) 

It is straightforward to deduce from Eq. (TT21 that the 
probability density function fr(t) is monotone decreasing 
from the level /r(0) = J °° n(x)p(dx) (which can be either 
finite or infinite) to the level /t(oo) = 0. 

Temporal hazard rate. The hazard rate of reaction 
time T is defined in terms of 



M*) 



lim — Pr(T < t + S\T > t) 

<5->0 



h{t) 



Pr(T > t) 



(13) 



with t > 0. Namely, the hazard rate hx(t) is the realiz- 
ation rate of the random variable T at time t, provided 
that it has not already occurred up to time t. The hazard 
rate is most commonly used in applied probability and in 
reliability theory |32h34| . Substituting the tail distribu- 
tion function Pr(T > t) of Eq. (ITU)) and the probability 
density function f T (t) of Eq. (JT2J into Eq. flTJJl yields the 
hazard rate of the reaction time T, 



h T (r) 



exp{-n(x)t}n(x)p(dx), t > 0. (14) 



It is evident from Eq. (|14p that the hazard rate 
hrit), analogous to the probability density function 
frit), is monotone decreasing from the level hr{0) — 
J °° rj{x)p{dx) (which can be either finite or infinite) to 
the level hr(oo) = 0. Note that a monotone decreas- 
ing hazard rate /ir(^) implies that the realization rate 
of the reaction time T diminishes in the course of time. 
Namely, the longer we are waiting for the reaction time 
T to realize, the lower the likelihood that it occurs im- 
mediately. 

Reaction range. The probability density function of 
the reaction range X is given by 

fx(x) = (l--E[exp{-ri(x)T}])//(x), x > 0. (15) 
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Eq. (|T5|) implies that if the Poissonian intensity p'{x) is 
a monotone decreasing function then so is the probab- 
ility density function fx(x). In general, however, the 
probability density function fx{x) itself is not necessar- 
ily monotone decreasing. Eq. f| 15[) further implies that 
the tail distribution function of the reaction range X is 
given by 

/>oo 

Pr(X >x)= (1 - E[exp{-r)(y)T}])p{dy), x > 0. 

J X 

(16) 

Note that both the probability density function fx{x) 
and the tail distribution function Pr(X > x) of the reac- 
tion range X involve the Laplace transform of the reac- 
tion time T. 

Reaction pair. The joint probability density function 
of the reaction pair (T, X) is given by 

f(T,x)(t,x) =Pr(T > t)exp{-tri(x)}ri{x)p'(x), t,x>0. 

(17) 

Comparing the joint probability density function 
f{T,x)(t, x) from Eq. (fTT|) with the product fT{t)fx(x) of 
the probability density function fr{t) from Eq. (IT2|) and 
the probability density function fx{x) from Eq. (|15|) . it 
is straightforward to see that 



f(T,x)(t,x)^f T (t)fx{x), t,x>0. 



(18) 



Relation (fT8[) implies that the reaction time T and the 
reaction range X are dependent random variables. This 
dependence is diametric to the conditional independence 
of the reaction time T and the reaction range X, con- 
ditioned on the realizations of distances {Di}, which is 
manifested by Eq. Thus the random Poissonian 

structure of the distances {Di} induces a statistical de- 
pendence between the reaction time T and the reaction 
range X. 

Conditional distributions. The conditional distri- 
bution of the reaction time T, conditioned on the realiz- 
ation of the reaction range X, is given by the tail distri- 
bution function 



Pr(T > t\X = x) 



j°° p r (y > s ) exp{-n(x)s}ds 
Jq 00 Pr(T > s) exp{— i](x)s}ds ' 



(19) 

for t > 0. Also, the conditional distribution of the reac- 
tion range X , conditioned to the realization of the reac- 
tion time T, is given by the tail distribution function 



Pt(X > x\T = t) = 



JT ex p{-v ( ,v) t }v(y)p( d y) 

g nftv exp{-n(y)t}n(y)p(dy) 



(20) 
for x > 0. 

As noted at the end of section II. A, the RARE model 
can be viewed as a descendant of Blumen's unified 
donor-acceptor energy transfer model [28j and Thomas 
et al.'s donor-acceptor recombination model [HI- In 
Refs. [Ill [5^1 the reaction time T was analyzed, and a 
counterpart of Eq. (ITUl) was established (in the context of 



Bernoulli scattering of excitations over an N-site lattice) . 
The analysis carried out in this section goes significantly 
beyond those of Refs. [H, [29[ as it addresses the mar- 
ginal, joint, and conditional distributions of the reaction 
pair (T, X). Moreover, all results established in this sec- 
tion are valid in the context of general Poisson scattering 
of excitations over general metric spaces, a setting which 
allows for high variability and robustness, whilst yielding 
closed-form analytic expressions. 



IV. MONTE-CARLO SIMULATION 

The stochastic analysis presented in the previous sec- 
tion facilitates the construction of a numerical Monte- 
Carlo algorithm for the simulation of the reaction pair 
(T,X). The steps of the Monte-Carlo algorithm are as 
follows: 

(1) Identify the inputs of the RARE model: the scat- 
tering function p(x) and the reactivity function rj(x). 

(2) Using the scattering function p(x) and the reactiv- 
ity function rj(x), numerically compute the integral 



I{t,x) 



e X p{-r,(y)t}r,{y)p(dy), t,x>0. (21) 



(3) Using the integral I(t, x), numerically compute the 
cumulative distribution function of the reaction time 

F(t) = 1 - exp (- I(t, oo)dT^j , t > 0, (22) 

and then numerically compute its inverse function 
F-^u) (0< u< 1). 

(4) Using the integral I(t, x), numerically compute the 
conditional cumulative distribution function of the reac- 
tion range 



G(x;t) = ^\, x>0, 
1 (r, oo) 



(23) 



and then numerically compute its inverse function 
G _1 (u; t) (0 < u < 1). 

(5) Generate a pair (1/1,1/2) of independent random 
variables, which are uniformly distributed over the unit 
interval, and generate the reaction pair (T,X) via 



T=F- 1 (Ui) and X = G _1 (C/ 2 ; 



(24) 



The third step of the Monte-Carlo algorithm is based 
on Eq. (flOl) . the fourth step is based on Eq. ([20]) . and in 
the fifth step we applied the following basic simulation 
principle Hi: If $(x) (x > 0) is the cumulative distribu- 
tion function of a positive- valued random variable £, and 
if U is a random variable which is uniformly distributed 
over the unit interval (and independent of £), then the 
random variable $ _1 (C/) is equal in law to the random 
variable £. 

The Monte-Carlo algorithm can be further used to sim- 
ulate random walks whose dynamics are governed by the 
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RARE model. To that end assume that the scattering of 
the excitations is spatially homogeneous. This homogen- 
eity implies that the Poissonian structure of the distances 
{Di} is invariant with respect to the position P of the re- 
action center. Now consider the following reaction-based 
propagation scheme: if the reaction center reacted with 
excitation i then the center jumps from its initial posi- 
tion P to the position of excitation i, Pi, and thereafter 
the process starts anew. This propagation scheme gener- 
ates a renewal Continuous Time Random Walk (CTRW) 
36 38], whose law of motion is as follows: 

(1) Initiate from an arbitrary position Pq at time to, 
and simulate a reaction pair (T, X). 

(2) At time t\ = to + T move to a point Pi which is 
uniformly distributed on a sphere with radius X centered 
at the point Po- 

(3) Set to := t\ and Po : = Pi, an d go back to step (1). 
The reaction time T is the CTRW's generic waiting 

time, and the reaction range X is CTRW's generic jump 
length. We emphasize that since the reaction time T and 
the reaction range X are dependent random variables, 
the RARE model induces CTRWs with coupled waiting 
times and jump lengths [39l - l4l| . 



V. POWER-LAW INPUTS 

To illustrate the results established so far we consider 
now the example of a RARE process with power-law in- 
puts. Specifically, we consider both the scattering func- 
tion and the reactivity function to be power-laws 



p{x) — ax a and i](x) — bx 13 , 



x>0 



(25) 



whose coefficients a and b as well as the exponents a and 
(3 are all positive parameters. 

Power-law inputs occur naturally in many cases. If 
we scatter the excitations uniformly across d-dimensional 
Euclidean space using a spatially homogeneous Poisson 
process, then the scattering function p(x) will be a power- 
law with exponent a — d. Indeed, if the excitations are 
scattered uniformly across the d-dimensional Euclidean 
space, then the mean number of excitations present in 
a ball of radius x, quantified by the scattering function 
p(x), is proportional to x d . Moreover, in many fractal 
settings, which are abundant in the physical sciences, the 
scattering of the chemical agents is uniform across some 
fractal unbounded subset of the d-dimensional Euclidean 
space. In such fractal settings the scattering function 
p(x) will be a power-law with exponent a, that equals 
the fractal dimension of the underlying fractal subset 
42]. On the other hand, power-law decays of physical 
and chemical interactions as a function of the distance 
between the interacting elements are prevalent in the the 
physical sciences, the best known examples being New- 
ton's law of gravitation, Coulomb's law, or the van der 
Waals law. 

We turn now to describe the statistical behavior of 
the RARE model with power-law inputs. For the RARE 



model to be well-defined, the exponent a of the scatter- 
ing function must be smaller than the exponent j3 of the 
reactivity function, and thus we assume that a < /3. 

Reaction time. From Eq. (TT0)) we obtain that the 
tail distribution function of the reaction time T is given 
by 



Pr(P >t) = exp [-cit a/f} ^ , t > 0, 



(26) 



where the precise value of the coefficient c\ is given by 
ci = r(l — a/ f3)ab a ^ . The tail distribution function 
of Eq. (|26|) characterizes a stretched exponential law (43l — 
|45| . The hazard rate of the reaction time T is a power-law 
with a negative exponent equal to —(1 — a/f3), and the 
moments of the reaction time T are given by 



E[T m ] = T(l + m/3/a)c 



i/3/a 



m > 0. 



(27) 



We emphasize that although the reaction time T has con- 
vergent moments of all orders, its moment generating 
function is divergent: E[exp(6T)] = oo for all 8 > 0. This 
statistical behavior of the stretched exponential distribu- 
tion, convergent moments and divergent moment gener- 
ating function, implies that the reaction time T displays 
a form of randomness which Mandelbrot categorized as 
'borderline randomness' [46|. Another well-known prob- 
ability distribution displaying such borderline random- 
ness is the log- normal distribution {ttT ]. 

Reaction range. Evaluating Eq. (fTB)) . while applying 
a moment-expansion of the Laplace transform of the reac- 
tion time T, we obtain that the tail distribution function 
of the reaction range X is given by the power-expansion 

Pr(A > x) = aa V (-l) m+1 \. H * \ , , 
^ x ' m\(mf3 - a) x m P- a 



rn—1 



(28) 

with x > 0. In turn, Eq. (|2"8")l implies that the asymptotic 
behavior of the tail distribution function of the reaction 
range X is given by 



Pr(X > x) 



C-2 



r fl- a 



, x — > oo, 



(29) 



where the coefficient C2 is given by C2 = [r(l + 
^/a)a 1 "^/ Q ]/[r(l-a//3)(/?-a)]. The tail asymptotics of 
Eq. (|29p characterize an as ymp totically Paretian distri- 
bution with exponent (3 — a [30, 48-5(J. The moments of 
the reaction range X are convergent, i.e., E[A m ] < oo, if 
and only if the exponent m is in the range < m < (3 — a. 
In particular, the reaction range X has a convergent 
mean if and only if 1 + a < (3. This statistical behavior 
of asymptotically Paretian distributions, convergent mo- 
ments only up to a given order, implies that the the reac- 
tion range X displays a form of randomness which Man- 
delbrot categorized as 'wild randomness' 0, [5l|. For 
a recent treatment of the 'categorization of randomness' 
see Ref. [H. 

Conditional distribution. From Eq. (|2"D|) we obtain 
that the conditional distribution of the reaction range X , 
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conditioned on the realization of the reaction time T, is 
given by the tail distribution function 



function 



Pr(A > x\T = t) = 



tbx 



exp(—u)u a / ^ 
T(l - a//3) 



du, x > 0. 

(30) 



Note that the integrand on the right hand side of Eq. (|30|) 
is the probability density function of a Gamma distri- 
bution with exponent 1 — a/(3. Thus, if we set £ to 
be a Gamma-distributed random variable with exponent 
1 - a/P, we obtain that Pr(X > x\T = t) = Pr(£ < 
tbx -13 ). Eq. ([50)1 further implies that the asymptotic be- 
havior of the conditional tail distribution function of the 
reaction range X is given by 



Pr(X >x\T = t)^c 3 



t i-a/p 



X — > OO, 



(31) 



where c 3 = 6 1 " Q /' 9 /r(2-a//3). As in the case of Eq. ([25]). 
the tail asymptotics of Eq. ([3~T1) characterize an asymp- 
totically Paretian distribution with exponent /3 — a. Con- 
sequently, also with any given realization of the reaction 
time T, the reaction range X displays wild randomness. 



VI. THERMODYNAMIC LIMIT 

In this section we explore the statistical behavior of 
the RARE model as the concentration of the excitations 
is increased to infinity. Specifically, we increase by n- 
fold the concentration of the excitations, and examine 
the limiting statistical behavior of the reaction time and 
the reaction range in the limit n — > oo. 

An n- fold increase of the concentration of the excita- 
tions results in replacing the 'original' scattering function 
p{x) by the 'concentrated' scattering function p n (x) = 
np(x). On the other hand, an ri-fold increase of the con- 
centration of the excitations does not affect the reactiv- 
ity function r/(x). In this section we denote by (T n ,X n ) 
the reaction pair corresponding to a RARE model with 
scattering function p n (x) and reactivity function rj(x). 
Moreover, we set 



rj(x)p(dx) 



(32) 



and assume that the integral on the right-hand-side of 
Eq. ([32]) is convergent. 

We introduce the scaled reaction time S n = nT n . A 
stochastic limit analysis detailed in Appendix E asserts 
that the joint probability density function of the reaction 
pair (S n ,X n ) attains the limit 

lim f(s n ,x n )(.s,x) = exp(-Xs)ri(x)p'(x), s, x > 0. 

Namely, the reaction pair (S n ,X n ) converges in law (as 
n — s- oo) to a limiting random vector (Soo,Xoo), whose 
distribution is governed by the joint probability density 



/(Soo,x 00 )(s,a ; ) = Xexp(-Xs)-r}(x)p'(x), s,x > 0. 

(34) 

Eq. implies the following: 

(1) The scaled reaction time S n converges in law to 
a stochastic limit Soo, which is exponentially distributed 
with rate A and density function 



fs^(s) = Aexp(-As), s > 0. 



(35) 



(2) The reaction range X n converges in law to a 
stochastic limit Xoo whose distribution is governed by 
the probability density function 



fXco( x ) = tV( x )p'( x )i x > 0. 



(36) 



(3) The joint probability density function of the ran- 
dom vector (Soo, Xoo) equals the product of its corres- 
ponding marginal probability density functions 



/(Soo.Xx,)^) = fs ao (s)fx oa {x), s,x> 0, 



(37) 



and hence the stochastic limits S 1 ,^ and Xoo are inde- 
pendent random variables. 

The stochastic limit result established in this section 
gets us 'all around the circle'. Indeed, in the preliminary 
analysis of section |TT] we obtained that the conditional 
distribution of the reaction time T, conditioned on the 
realizations of distances {-Dj}, is exponential and inde- 
pendent of the reaction range X. However, in section ITTT1 
we found that the distribution of the reaction time T is 
general, and is tightly coupled to the distribution of the 
reaction range X. Thus, the random Poissonian struc- 
ture of the distances {Di} shifts the distribution of the 
reaction time T from exponential to general, and induces 
a dependence between the reaction time T and the reac- 
tion range X. In this section we established that in the 
infinite concentration limit n — > oo the original statistical 
structure is recovered: (i) the exponential distribution is 
restored, as the stochastic limit of the scaled reaction 
time is exponentially distributed; and (ii) the independ- 
ence is restored, as the stochastic limit Soo 01 the scaled 
reaction time and the stochastic limit X^ of the reaction 
range are independent random variables. 

To illustrate the stochastic limit result established in 
this section consider the example of power law scattering 
functions studied in the previous section and exponential 
reactivity functions 

p(x) = ax a and rj(x) = bexp(—(3x), x > 0, (38) 

where the coefficients a and b as well as the parameters 
a and /3, are all positive. In this example the stochastic 
limit Soo is exponentially distributed with rate 



A = ah 



r(l + q) 

f3 a ' 



(39) 



and the stochastic limit is Gamma distributed with 
probability density function 



fit* 



exp(-f3x)x a , x > 



(40) 



Note that the Gamma probability density function of 
Eq. (|40l) is (i) unbounded and monotone decreasing in the 
exponent range a < 1; (ii) bounded and monotone de- 
creasing at the exponent value a = 1, in which case Xoo 
is exponentially distributed with mean 1//3; (hi) bounded 
and unimodal, with mode x* = (a — l)//3, in the expo- 
nent range a > 1. 



VII. CONCLUSION 

In this paper we presented a general spatio-chemical 
stochastic model for generalized RAndom Relaxation 
(RARE) processes in complex disordered systems. The 
RARE model considers a collection of excitations which 
are randomly scattered around a reaction center in some 
general embedding metric space. The RARE model has 
two input quantities: (i) the scattering function p(x) 
quantifying the scattering intensity of the excitations 
around the reaction center as a function of the distance 
(x > 0) from the center; (ii) the reactivity function rj{x) 
quantifying the reaction rate between the excitations and 
the reaction center as a function of the distance (x > 0) 
from the center. The scattering and reactivity functions 
provide a straightforward intuitive description, as well as 
a precise mathematical formulation, of general relaxation 
processes in complex disordered systems. 

The RARE model has two random outputs: (i) the re- 
action time T of its random relaxation process; (ii) the 
reaction range X of its random relaxation process. A 
detailed stochastic analysis of the reaction pair (T, X) 
was carried out, yielding closed form results regarding 
the statistics of this pair: marginal distributions, joint 
distribution, and conditional distributions. The results 
established further led to a Monte-Carlo algorithm for 
the simulation of the model's random relaxation process, 
and to the thermodynamic of the RARE model. In ad- 
dition, we investigated in detail the case of power-law 
inputs, which were shown to yield stretched exponential 
relaxation patterns and asymptotically Paretian relaxa- 
tion ranges. 

The RARE model is a compact and transparent 
stochastic approach to non-exponential relaxation pro- 
cesses. On the one hand the model's inputs are both 
intuitively clear and mathematically precise, and are dir- 
ectly related to the physical properties of the system con- 
sidered. On the other hand the model is rather robust 
and versatile, as its general mathematical formulation 
accommodates diverse physical settings. We therefore 
expect the RARE model to be useful to a wide range of 
applications in the physical sciences and beyond. 
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Appendix A: The distribution of the reaction time T 

Eq. ([1} implies that the tail distribution function of 
the reaction time T is given by 

Pr(T > t) = Pr (mm{Ti} > tj = Pr (f){Ti > t}j . 

(Al) 

After conditioning with respect to the distances {A}, 
this is equal to 



Pr(T > t) = E 



Pr ( f|{T, > *}|{A} 



(A2) 



Using the assumptions of the RARE model and Eq. ((3])), 
we have 



Pr(T > t) = E 



E 



Y[PT{Ti>t\Di) 

, i 

l[vq>{-V(Di)t} 



(A3) 



With Campbell's theorem of the theory of Poisson pro- 
cesses (see Section 3.2 in Ref. [24j). 

Pr(T >t) = cxp f-j [1 - exp{-T)(x)t}]p(dx) 



Eqs. (|A"Tj) to ([AH) prove Eq. (ITUl) . 



(A4) 



Appendix B: The distribution of the reaction 
pair (T,X) 

In what follows we set 

E(t,9,x) = E[exp{-tR(0)-dR(x)}], t,9,x>0. (Bl) 
Eq. ([S]) implies that 

tR(0) + 9R(x) = {trj{Di)l{D t < x) 

i 

+(i + %(A)I(A >xj), (B2) 
and hence Eq. (|B1[) further implies that 

E(t,9,x) = E[jJexp(-try(A)I(A <x) 

i 

-{t + 0)ri(D i )l(D i >x))\ (B3) 
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Consequently, Campbell's theorem of the theory of Pois- where 
son processes (see Section 3.2 in Ref. (24|) implies that 

E(t, 9, x) = exp{-*(i, 9, x)}, (B4) 



1 - exp{-tr](s)I(s < x) - (t + 9)rj(s)I(s > x)} p(ds) 
1 - exp{-trj(s)} p(ds)+ / 1 -exp{-(t + 6)ri(s)} p{ds). 

J X 

I 



(B5) 



Applying conditioning with respect to the distances 
{Di} and using Eq. @, we obtain that the joint tail 
distribution function of the reaction pair (T, X) is given 
by 

Pr(T > t,X > x) = E[Pr(T > t,X > x\{D z })} 

R{x) 



E 



exp{-ti?(0)} 



R(0) 



Differentiating this expression with respect to the vari- 
able t yields 



Pr(T > t, X > x) = -~E[exp{-tR(0)}R(x)]. (B7) 



d_ 

m 



Conversely, differentiating Eq. (fBTj) with respect to the 
variable 9 yields 

^E{t, 8, x) = -E[ex P {-ti?(0) - 9R(x)}R(x)}. (B8) 
o9 

Combination of Eqs. (|B7| and (|B8|) yields 



-Px(T>t,X>x)= ^E(t,9,x) 



8=0 



Now, Eq. (|B4|) implies that 



)=0 



= -ex P {-*(i,0,x)} —m{t,9,x) 



(BIO) 



t, x > 0. (B6) and Eq. (TB5|) further implies that 



*(t,0,x) 



1 - exp{-*77(s)} p(ds), (Bll) 



and 



09 



V(t,9,x) 



>=0 J x 



3xp{-tr)(s)}r)(s)p(ds). (B12) 



(B9) Thus, combining Eqs. (|B9]) to (lB~T2)) . we conclude that 



- — Pr(T >t,X > x) = exp I- J \l - exp{-tr)(s)}\p(ds)J J exp{-trj(s)}r)(s)p(ds). (B13) 

Finally, differentiating this last expression with respect to the variable x and using Eq. (| 10[) , we arrive at the desired 
result, 



if 



dtdx 

such that Eq. (|B14|) proves Eq. (fI7]l 



Pr(T > t,X > x) = Pr(T > t) exp{-tr](x)}rj(x)p'(x), 



(B14) 



r 



Appendix C: The marginal distributions of the ity function of the random pair (T, X) over the distance 

reaction pair (T, X) variable x, 



The probability density function of the reaction time 
T is attained by integrating the joint probability dens- 



h(t) = / f {T .x)(t,x)dx. 
Jo 



(Cf) 
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With Eq. (HTJ), 

f T (t) =Pr(T >t) / cxp{-t?7(a;)}r/(a;)p(dx), t > 0. 
Jo 

(C2) 

Eqs. JCl]) and ([C2j prove Eq. ([T2]>. Note that the Laplace 
transform of the reaction time T is given by 

/>oo 

E[exp(-0T)] = / cxp{-6t)f T (t)dt. (C3) 
Jo 

Integration by parts yields 

/>oo 

E[exp(-0T)] =1-8 exp(-6t) Pr(T > (C4) 
Jo 

where 9 > 0. 

The probability density function of the reaction range 
X is attained by integrating the joint probability density 
function of the random pair (T, X) over the time variable 
t, 



(C5) 



= / f(r,x)(t,x)dt 



and, using Eq. (fT7|) . 

= Pr(T > t) e>q>{-t»j(z)}d^ 77(z)p' (*)■ 

(C6) 

With Eqs. JC3| and (CH) and (9 = we find 

/x(a;) = (l-E[exp{-7 ? (x)T}])p'(a ; ) ! x > 0. (C7) 
Eqs. ([C5]) to fCT} prove Eq. fig]) . 



Appendix D: The conditional distributions of the 
reaction pair (T, X) 

The conditional probability density function of the re- 
action time T, conditioned to the realization of the reac- 
tion range X, is given by 



fc\x=x(t) 



f(T,X)(t,x) 

fx(x) 



(Dl) 



With Eqs. (fT7|) and (fl5)). 

Pr(T > t) exp{-tri(x)}r)(x)p'(x) 



fr\x=x(t) 



(Io° P r (^ > T ) ex P{~ Tr](x)}dT)r](x)p'(x) 
Pr(T > t) exp{-ti](x)} 



r o °° Pr(T > r) exp{-rr](x)}dT ' 



(D2) 



In turn, Eqs. (|DT|) and (|D2I) imply Eq. (fT9|) . 

The conditional probability density function of the re- 
action range X, conditioned to the realization of the re- 
action time T, is given by 



fx\T=t{x) 



f(T,X)(t,x) 

Mt) ' 



(D3) 



Using Eqs. (JT7|) and ([T2]). we see that 

Pr(T > t)exp{-tr)(x)}r)(x)p'(x) 



fx\T=t{x) = 



Pr(T > t) exp{-ri(x)t}ri(x)p(dx) 
exp{— tr/(x)}i](x)p' (x) 



JiT exp{-?7(a;)t}?7(x)p(da;) 
In turn, Eqs. |D3| and (JD4J) imply Eq. 



(D4) 



Appendix E: The limit distribution of the reaction 
pair (5„,X n ) 

The scaling S n = nT n implies the following connec- 
tion between the joint cumulative distribution function 
of the reaction pair (S n ,X n ) and the joint cumulative 
distribution function of the reaction pair (T n ,X n ), 



Pr(5„ < s, X n <x)= Pr (T„ < -, X n < xj , s, a; > 0. 

(El) 

Consequently, differentiating this expression with respect 
to the variables s and x yields the following connection 
between the joint probability density function of the reac- 
tion pair (S n , X n ) and the joint probability density func- 
tion of the reaction pair (T n ,X n ), 



f(S n ,x n )(s,x) = f(r n ,x n ) (-, x ) -■ 

Vn In 



(E2) 



Using Eq. (fTT]) with the concentrated scattering function 
p n [x)t we arrive at 



f(s n ,x n ){s,x) = Pr (r n > exp (--r)(x)J 



v(x)p' n (x) 
n 

(E3) 

With the definition of the concentrated scattering func- 
tion p„(x), this implies 

f(S n ,x n )(s,x) = Pr (T n > -) exp f--7?(x)J n(x)p'(x). 

(E4) 

Eqs. (1E2]) to ([E4]) lead us to conclude that 



lim Pr (T n > - 

n— >oo \ n 



lim fis n ,x n ){s,x) 



n(x)p'(x), 
(E5) 

with s, x > 0. Now, using Eq. (| 10[) with the concentrated 
scattering function p n {x) we find 



lim Pr 

n— >oo 



1 - exp (-rj(x)-) Pn(^) J ■ (E6) 



= lim exp 



With the definition of the concentrated scattering func- 
tion p n {x) this leads to 



lim Pr 

= exp 



lim 



1 — exp (-n(x)t /n) 
1/n 



p(dx) . (E7) 
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Finally, with de l'Hospital's rule and Eq. (|32[) we obtain Combining together the limits of Eqs. (IE4I) and (IE8[) 

proves Eq. 



lim Pr 

n— >ao \ n 



( r ">j) 

exp(- / [r)(x)i\p(dx) ) = exp(-At). (E8) 
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